# knitr::opts_chunk$set(echo = TRUE)
# knitr::opts_chunk$set(render = 'normal_print')
# knitr::opts_chunk$set(results='asis')

library(table1)
## 
## Attaching package: 'table1'
## The following objects are masked from 'package:base':
## 
##     units, units<-
library(ggplot2)
library(ggpubr)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggsci)
library(jtools)
library(lmerTest)
## Loading required package: lme4
## Loading required package: Matrix
## Registered S3 methods overwritten by 'lme4':
##   method                          from
##   cooks.distance.influence.merMod car 
##   influence.merMod                car 
##   dfbeta.influence.merMod         car 
##   dfbetas.influence.merMod        car
## 
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
## 
##     lmer
## The following object is masked from 'package:stats':
## 
##     step
df <- read.table(file = file.path(getwd(),"R_Input_2.txt"), 
                 header = TRUE, 
                 quote = "", 
                 sep = "\t", 
                 row.names = 6)
#my_comparisons <- list( c("NORMAL", "PMS2"), c("NORMAL", "MSH6"), c("MSH6", "MSH2_MSH6") )

###################################################################
# plot genomic landscape by MMRd types
####################################################################
# add in cancer types
df_Tumor.purity <- read.table(file = file.path(getwd(),"Tumor_purity.txt"), sep = "\t", header = T, quote = "" )
row.names(df_Tumor.purity) <- df_Tumor.purity$Sample_ID

df_Cancer <- read.table(file = "Cancer_type.txt", header = T, sep = "\t")
row.names(df_Cancer) <- df_Cancer$Sample_ID

df_1 <- merge(df, df_Cancer, by = 0)
row.names(df_1) <- df_1$Row.names
df_1 <- df_1[,-which(names(df_1) %in% c("Row.names","Sample_ID"))]

df_1 <- merge(df_1, df_Tumor.purity, by = 0)
row.names(df_1) <- df_1$Row.names
df_1 <- df_1[,-which(names(df_1) %in% c("Row.names","Sample_ID"))]

# Add the column "Multiplex_or_not"
df_2 <- df_1 %>% mutate(Multiplex_or_not = 
                  ifelse(df_1$CLASSIFICATION_PLOT %in% c("MSH2_MSH6","MLH1_PMS2","Complex"), "Multiplex",
                         ifelse(df_1$CLASSIFICATION_PLOT %in% c("MSH6","PMS2"), "Singleton",df_1$CLASSIFICATION_PLOT)))

# Add the column "MSH6_or_not"
df_2 <- df_2 %>% mutate(MSH6_or_not = 
                           ifelse(grepl("MSH6",df_2$CLASSIFICATION_PLOT),"MSH6","Not_MSH6"))

# Add the column "PMS2_or_not"
df_2 <- df_2 %>% mutate(PMS2_or_not = 
                           ifelse(grepl("PMS2",df_2$CLASSIFICATION_PLOT),"PMS2","Not_PMS2"))

df_2 <- df_2 %>% mutate(Multiplex_or_not = as.factor(Multiplex_or_not)) %>% mutate(Multiplex_or_not = relevel(Multiplex_or_not, ref = "Singleton"))
df_2 <- df_2 %>% mutate(MSH6_or_not = as.factor(MSH6_or_not)) %>% mutate(MSH6_or_not = relevel(MSH6_or_not, ref = "MSH6"))
df_2 <- df_2 %>% mutate(PMS2_or_not = as.factor(PMS2_or_not)) %>% mutate(PMS2_or_not = relevel(PMS2_or_not, ref = "PMS2"))
df_2 <- df_2 %>% mutate(Cancer_Type_NEW = as.factor(Cancer_Type_NEW)) %>% mutate(Cancer_Type_NEW = relevel(Cancer_Type_NEW, ref = "Colon"))

row.names(df_2) <- row.names(df_1)

# renames CLASSIFICATION_PLOT to MMRd_Type for both df_1 and df_2
df_1 <- df_1 %>% rename(MMRd_Type = CLASSIFICATION_PLOT) 
df_2 <- df_2 %>% rename(MMRd_Type = CLASSIFICATION_PLOT) 

label(df_1$Current_Age) <- "Age"
units(df_1$Current_Age) <- "years"

label(df_2$Current_Age) <- "Age"
units(df_2$Current_Age) <- "years"

label(df_1$Tumor_Purity) <- "Tumor Purity"
units(df_1$Tumor_Purity) <- "%"

label(df_2$Tumor_Purity) <- "Tumor Purity"
units(df_2$Tumor_Purity) <- "%"

label(df_1$Cancer_Type_NEW) <- "Cancer Type"
label(df_2$Cancer_Type_NEW) <- "Cancer Type"

label(df_1$Sample_Type) <- "Sample Type"
label(df_2$Sample_Type) <- "Sample Type"
###################################################################
# make table1 for both df_1 and df_2 for baseline characteristic table
###################################################################
#topclass="Rtable1-grid Rtable1-shade Rtable1-times"
table1(~ Current_Age + Race + Ethnicity + Tumor_Purity + Cancer_Type_NEW + Sample_Type | MMRd_Type, data=df_1[df_1$MMRd_Type!="Delete",], overall="Total")
?
(N=5)		 Complex
(N=11)		 MLH1_PMS2
(N=251)		 MSH2_MSH6
(N=70)		 MSH6
(N=12)		 NORMAL
(N=19)		 PMS2
(N=13)		 Total
(N=381)
Current_Age
Mean (SD) 64.4 (12.8) 63.8 (19.9) 67.9 (13.2) 56.8 (14.1) 64.5 (11.2) 61.1 (14.5) 54.3 (16.7) 64.8 (14.4)
Median [Min, Max] 63.0 [46.0, 77.0] 70.0 [26.0, 88.0] 71.0 [25.0, 90.0] 55.0 [30.0, 90.0] 67.0 [36.0, 76.0] 65.0 [31.0, 87.0] 48.0 [37.0, 84.0] 67.0 [25.0, 90.0]
Race
BLACK OR AFRICAN AMERICAN 1 (20.0%) 0 (0%) 15 (6.0%) 4 (5.7%) 0 (0%) 1 (5.3%) 1 (7.7%) 22 (5.8%)
PT REFUSED TO ANSWER 1 (20.0%) 0 (0%) 6 (2.4%) 1 (1.4%) 0 (0%) 1 (5.3%) 0 (0%) 9 (2.4%)
WHITE 3 (60.0%) 8 (72.7%) 208 (82.9%) 57 (81.4%) 10 (83.3%) 16 (84.2%) 9 (69.2%) 311 (81.6%)
ASIAN-FAR EAST/INDIAN SUBCONT 0 (0%) 2 (18.2%) 11 (4.4%) 4 (5.7%) 1 (8.3%) 1 (5.3%) 3 (23.1%) 22 (5.8%)
OTHER 0 (0%) 1 (9.1%) 9 (3.6%) 4 (5.7%) 1 (8.3%) 0 (0%) 0 (0%) 15 (3.9%)
NATIVE AMERICAN-AM IND/ALASKA 0 (0%) 0 (0%) 2 (0.8%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 2 (0.5%)
Ethnicity
Non-Spanish; Non-Hispanic 5 (100%) 8 (72.7%) 215 (85.7%) 58 (82.9%) 11 (91.7%) 18 (94.7%) 11 (84.6%) 326 (85.6%)
"Spanish NOS; Hispanic NOS, Latino NOS" 0 (0%) 2 (18.2%) 14 (5.6%) 8 (11.4%) 1 (8.3%) 1 (5.3%) 1 (7.7%) 27 (7.1%)
Dominican Republic 0 (0%) 1 (9.1%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 1 (0.3%)
0 0 (0%) 0 (0%) 3 (1.2%) 2 (2.9%) 0 (0%) 0 (0%) 0 (0%) 5 (1.3%)
Puerto Rican 0 (0%) 0 (0%) 2 (0.8%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 2 (0.5%)
Unknown whether Spanish or not 0 (0%) 0 (0%) 17 (6.8%) 2 (2.9%) 0 (0%) 0 (0%) 1 (7.7%) 20 (5.2%)
Tumor_Purity
Mean (SD) 48.0 (19.2) 38.2 (12.7) 38.3 (16.0) 39.6 (15.8) 39.6 (16.3) 47.2 (17.0) 38.8 (19.4) 39.2 (16.1)
Median [Min, Max] 50.0 [20.0, 70.0] 30.0 [25.0, 60.0] 40.0 [10.0, 90.0] 40.0 [5.00, 80.0] 30.0 [25.0, 80.0] 50.0 [20.0, 70.0] 35.0 [10.0, 70.0] 40.0 [5.00, 90.0]
Missing 0 (0%) 0 (0%) 1 (0.4%) 2 (2.9%) 0 (0%) 0 (0%) 0 (0%) 3 (0.8%)
Cancer_Type_NEW
Bladder 1 (20.0%) 0 (0%) 0 (0%) 5 (7.1%) 2 (16.7%) 0 (0%) 1 (7.7%) 9 (2.4%)
Endometrial 1 (20.0%) 1 (9.1%) 62 (24.7%) 14 (20.0%) 1 (8.3%) 4 (21.1%) 0 (0%) 83 (21.8%)
Esophagogastric 2 (40.0%) 1 (9.1%) 33 (13.1%) 4 (5.7%) 0 (0%) 2 (10.5%) 0 (0%) 42 (11.0%)
Small Bowel 1 (20.0%) 0 (0%) 4 (1.6%) 2 (2.9%) 0 (0%) 1 (5.3%) 4 (30.8%) 12 (3.1%)
Colon 0 (0%) 9 (81.8%) 118 (47.0%) 19 (27.1%) 4 (33.3%) 7 (36.8%) 4 (30.8%) 161 (42.3%)
Biliary 0 (0%) 0 (0%) 5 (2.0%) 0 (0%) 0 (0%) 0 (0%) 1 (7.7%) 6 (1.6%)
Colorectal 0 (0%) 0 (0%) 10 (4.0%) 5 (7.1%) 3 (25.0%) 4 (21.1%) 1 (7.7%) 23 (6.0%)
CUP 0 (0%) 0 (0%) 5 (2.0%) 0 (0%) 0 (0%) 1 (5.3%) 0 (0%) 6 (1.6%)
Other 0 (0%) 0 (0%) 2 (0.8%) 5 (7.1%) 0 (0%) 0 (0%) 0 (0%) 7 (1.8%)
Ovarian 0 (0%) 0 (0%) 1 (0.4%) 2 (2.9%) 0 (0%) 0 (0%) 1 (7.7%) 4 (1.0%)
Pancreas 0 (0%) 0 (0%) 1 (0.4%) 2 (2.9%) 1 (8.3%) 0 (0%) 0 (0%) 4 (1.0%)
Prostate 0 (0%) 0 (0%) 1 (0.4%) 2 (2.9%) 1 (8.3%) 0 (0%) 0 (0%) 4 (1.0%)
Rectal 0 (0%) 0 (0%) 9 (3.6%) 10 (14.3%) 0 (0%) 0 (0%) 1 (7.7%) 20 (5.2%)
Sample_Type
Metastasis 1 (20.0%) 1 (9.1%) 29 (11.6%) 6 (8.6%) 2 (16.7%) 4 (21.1%) 2 (15.4%) 45 (11.8%)
Mixed 1 (20.0%) 1 (9.1%) 17 (6.8%) 5 (7.1%) 2 (16.7%) 1 (5.3%) 1 (7.7%) 28 (7.3%)
Primary 3 (60.0%) 9 (81.8%) 205 (81.7%) 59 (84.3%) 8 (66.7%) 14 (73.7%) 10 (76.9%) 308 (80.8%) 
table1(~ Current_Age + Race + Ethnicity + Tumor_Purity + Cancer_Type_NEW + Sample_Type | Multiplex_or_not, data=df_2[df_2$MMRd_Type!="Delete",], overall="Total")
Singleton
(N=25)		 ?
(N=5)		 Multiplex
(N=332)		 NORMAL
(N=19)		 Total
(N=381)
Current_Age
Mean (SD) 59.2 (15.0) 64.4 (12.8) 65.4 (14.3) 61.1 (14.5) 64.8 (14.4)
Median [Min, Max] 64.0 [36.0, 84.0] 63.0 [46.0, 77.0] 68.0 [25.0, 90.0] 65.0 [31.0, 87.0] 67.0 [25.0, 90.0]
Race
ASIAN-FAR EAST/INDIAN SUBCONT 4 (16.0%) 0 (0%) 17 (5.1%) 1 (5.3%) 22 (5.8%)
BLACK OR AFRICAN AMERICAN 1 (4.0%) 1 (20.0%) 19 (5.7%) 1 (5.3%) 22 (5.8%)
OTHER 1 (4.0%) 0 (0%) 14 (4.2%) 0 (0%) 15 (3.9%)
WHITE 19 (76.0%) 3 (60.0%) 273 (82.2%) 16 (84.2%) 311 (81.6%)
PT REFUSED TO ANSWER 0 (0%) 1 (20.0%) 7 (2.1%) 1 (5.3%) 9 (2.4%)
NATIVE AMERICAN-AM IND/ALASKA 0 (0%) 0 (0%) 2 (0.6%) 0 (0%) 2 (0.5%)
Ethnicity
"Spanish NOS; Hispanic NOS, Latino NOS" 2 (8.0%) 0 (0%) 24 (7.2%) 1 (5.3%) 27 (7.1%)
Non-Spanish; Non-Hispanic 22 (88.0%) 5 (100%) 281 (84.6%) 18 (94.7%) 326 (85.6%)
Unknown whether Spanish or not 1 (4.0%) 0 (0%) 19 (5.7%) 0 (0%) 20 (5.2%)
0 0 (0%) 0 (0%) 5 (1.5%) 0 (0%) 5 (1.3%)
Dominican Republic 0 (0%) 0 (0%) 1 (0.3%) 0 (0%) 1 (0.3%)
Puerto Rican 0 (0%) 0 (0%) 2 (0.6%) 0 (0%) 2 (0.5%)
Tumor_Purity
Mean (SD) 39.2 (17.6) 48.0 (19.2) 38.6 (15.8) 47.2 (17.0) 39.2 (16.1)
Median [Min, Max] 30.0 [10.0, 80.0] 50.0 [20.0, 70.0] 40.0 [5.00, 90.0] 50.0 [20.0, 70.0] 40.0 [5.00, 90.0]
Missing 0 (0%) 0 (0%) 3 (0.9%) 0 (0%) 3 (0.8%)
Cancer_Type_NEW
Colon 8 (32.0%) 0 (0%) 146 (44.0%) 7 (36.8%) 161 (42.3%)
Biliary 1 (4.0%) 0 (0%) 5 (1.5%) 0 (0%) 6 (1.6%)
Bladder 3 (12.0%) 1 (20.0%) 5 (1.5%) 0 (0%) 9 (2.4%)
Colorectal 4 (16.0%) 0 (0%) 15 (4.5%) 4 (21.1%) 23 (6.0%)
CUP 0 (0%) 0 (0%) 5 (1.5%) 1 (5.3%) 6 (1.6%)
Endometrial 1 (4.0%) 1 (20.0%) 77 (23.2%) 4 (21.1%) 83 (21.8%)
Esophagogastric 0 (0%) 2 (40.0%) 38 (11.4%) 2 (10.5%) 42 (11.0%)
Other 0 (0%) 0 (0%) 7 (2.1%) 0 (0%) 7 (1.8%)
Ovarian 1 (4.0%) 0 (0%) 3 (0.9%) 0 (0%) 4 (1.0%)
Pancreas 1 (4.0%) 0 (0%) 3 (0.9%) 0 (0%) 4 (1.0%)
Prostate 1 (4.0%) 0 (0%) 3 (0.9%) 0 (0%) 4 (1.0%)
Rectal 1 (4.0%) 0 (0%) 19 (5.7%) 0 (0%) 20 (5.2%)
Small Bowel 4 (16.0%) 1 (20.0%) 6 (1.8%) 1 (5.3%) 12 (3.1%)
Sample_Type
Metastasis 4 (16.0%) 1 (20.0%) 36 (10.8%) 4 (21.1%) 45 (11.8%)
Mixed 3 (12.0%) 1 (20.0%) 23 (6.9%) 1 (5.3%) 28 (7.3%)
Primary 18 (72.0%) 3 (60.0%) 273 (82.2%) 14 (73.7%) 308 (80.8%) 
table1(~ Current_Age + Race + Ethnicity + Tumor_Purity + Cancer_Type_NEW + Sample_Type | MSH6_or_not, data=df_2[df_2$MMRd_Type!="Delete",], overall="Total")
MSH6
(N=82)		 Not_MSH6
(N=299)		 Total
(N=381)
Current_Age
Mean (SD) 57.9 (13.9) 66.6 (14.0) 64.8 (14.4)
Median [Min, Max] 56.5 [30.0, 90.0] 70.0 [25.0, 90.0] 67.0 [25.0, 90.0]
Race
ASIAN-FAR EAST/INDIAN SUBCONT 5 (6.1%) 17 (5.7%) 22 (5.8%)
BLACK OR AFRICAN AMERICAN 4 (4.9%) 18 (6.0%) 22 (5.8%)
OTHER 5 (6.1%) 10 (3.3%) 15 (3.9%)
PT REFUSED TO ANSWER 1 (1.2%) 8 (2.7%) 9 (2.4%)
WHITE 67 (81.7%) 244 (81.6%) 311 (81.6%)
NATIVE AMERICAN-AM IND/ALASKA 0 (0%) 2 (0.7%) 2 (0.5%)
Ethnicity
"Spanish NOS; Hispanic NOS, Latino NOS" 9 (11.0%) 18 (6.0%) 27 (7.1%)
0 2 (2.4%) 3 (1.0%) 5 (1.3%)
Non-Spanish; Non-Hispanic 69 (84.1%) 257 (86.0%) 326 (85.6%)
Unknown whether Spanish or not 2 (2.4%) 18 (6.0%) 20 (5.2%)
Dominican Republic 0 (0%) 1 (0.3%) 1 (0.3%)
Puerto Rican 0 (0%) 2 (0.7%) 2 (0.5%)
Tumor_Purity
Mean (SD) 39.6 (15.8) 39.1 (16.2) 39.2 (16.1)
Median [Min, Max] 40.0 [5.00, 80.0] 40.0 [10.0, 90.0] 40.0 [5.00, 90.0]
Missing 2 (2.4%) 1 (0.3%) 3 (0.8%)
Cancer_Type_NEW
Colon 23 (28.0%) 138 (46.2%) 161 (42.3%)
Biliary 0 (0%) 6 (2.0%) 6 (1.6%)
Bladder 7 (8.5%) 2 (0.7%) 9 (2.4%)
Colorectal 8 (9.8%) 15 (5.0%) 23 (6.0%)
CUP 0 (0%) 6 (2.0%) 6 (1.6%)
Endometrial 15 (18.3%) 68 (22.7%) 83 (21.8%)
Esophagogastric 4 (4.9%) 38 (12.7%) 42 (11.0%)
Other 5 (6.1%) 2 (0.7%) 7 (1.8%)
Ovarian 2 (2.4%) 2 (0.7%) 4 (1.0%)
Pancreas 3 (3.7%) 1 (0.3%) 4 (1.0%)
Prostate 3 (3.7%) 1 (0.3%) 4 (1.0%)
Rectal 10 (12.2%) 10 (3.3%) 20 (5.2%)
Small Bowel 2 (2.4%) 10 (3.3%) 12 (3.1%)
Sample_Type
Metastasis 8 (9.8%) 37 (12.4%) 45 (11.8%)
Mixed 7 (8.5%) 21 (7.0%) 28 (7.3%)
Primary 67 (81.7%) 241 (80.6%) 308 (80.8%) 
table1(~ Current_Age + Race + Ethnicity + Tumor_Purity + Cancer_Type_NEW + Sample_Type | PMS2_or_not, data=df_2[df_2$MMRd_Type!="Delete",], overall="Total")
PMS2
(N=264)		 Not_PMS2
(N=117)		 Total
(N=381)
Current_Age
Mean (SD) 67.2 (13.7) 59.3 (14.6) 64.8 (14.4)
Median [Min, Max] 70.0 [25.0, 90.0] 59.0 [26.0, 90.0] 67.0 [25.0, 90.0]
Race
ASIAN-FAR EAST/INDIAN SUBCONT 14 (5.3%) 8 (6.8%) 22 (5.8%)
BLACK OR AFRICAN AMERICAN 16 (6.1%) 6 (5.1%) 22 (5.8%)
NATIVE AMERICAN-AM IND/ALASKA 2 (0.8%) 0 (0%) 2 (0.5%)
OTHER 9 (3.4%) 6 (5.1%) 15 (3.9%)
PT REFUSED TO ANSWER 6 (2.3%) 3 (2.6%) 9 (2.4%)
WHITE 217 (82.2%) 94 (80.3%) 311 (81.6%)
Ethnicity
"Spanish NOS; Hispanic NOS, Latino NOS" 15 (5.7%) 12 (10.3%) 27 (7.1%)
0 3 (1.1%) 2 (1.7%) 5 (1.3%)
Non-Spanish; Non-Hispanic 226 (85.6%) 100 (85.5%) 326 (85.6%)
Puerto Rican 2 (0.8%) 0 (0%) 2 (0.5%)
Unknown whether Spanish or not 18 (6.8%) 2 (1.7%) 20 (5.2%)
Dominican Republic 0 (0%) 1 (0.9%) 1 (0.3%)
Tumor_Purity
Mean (SD) 38.3 (16.1) 41.1 (16.0) 39.2 (16.1)
Median [Min, Max] 40.0 [10.0, 90.0] 40.0 [5.00, 80.0] 40.0 [5.00, 90.0]
Missing 1 (0.4%) 2 (1.7%) 3 (0.8%)
Cancer_Type_NEW
Colon 122 (46.2%) 39 (33.3%) 161 (42.3%)
Biliary 6 (2.3%) 0 (0%) 6 (1.6%)
Bladder 1 (0.4%) 8 (6.8%) 9 (2.4%)
Colorectal 11 (4.2%) 12 (10.3%) 23 (6.0%)
CUP 5 (1.9%) 1 (0.9%) 6 (1.6%)
Endometrial 62 (23.5%) 21 (17.9%) 83 (21.8%)
Esophagogastric 33 (12.5%) 9 (7.7%) 42 (11.0%)
Other 2 (0.8%) 5 (4.3%) 7 (1.8%)
Ovarian 2 (0.8%) 2 (1.7%) 4 (1.0%)
Pancreas 1 (0.4%) 3 (2.6%) 4 (1.0%)
Prostate 1 (0.4%) 3 (2.6%) 4 (1.0%)
Rectal 10 (3.8%) 10 (8.5%) 20 (5.2%)
Small Bowel 8 (3.0%) 4 (3.4%) 12 (3.1%)
Sample_Type
Metastasis 31 (11.7%) 14 (12.0%) 45 (11.8%)
Mixed 18 (6.8%) 10 (8.5%) 28 (7.3%)
Primary 215 (81.4%) 93 (79.5%) 308 (80.8%) 
knitr::knit_exit()

plotCorr.MMRd <- function(x, y, z) {
  my_comparisons <- list()
  df_corr <- compare_means(as.formula(paste0(y," ~ ", z)),  data = x[-which(x$MMRd_Type=="Delete"),])
  idx <- 1
  for(sig in which(df_corr$p.signif != "ns")){
    #print(as.character(df_corr[sig,c("group1")]))
    
    my_comparisons[[idx]] <- c(as.character(df_corr[sig,c("group1")]), as.character(df_corr[sig,c("group2")]))
    idx <- idx + 1
  }
  n <- length(my_comparisons)
  write.table(df_corr, 
              file = file.path(getwd(),paste0(y,".",z,".corr.result.txt")), 
              quote = FALSE, 
              sep = "\t")
  
  ylabel <- y
  
  if (y == "Fraction_Genome_Altered") {
    max_y <- max(x[,y])
    min_y <- min(x[,y])
    assign(paste0("p_",y), ggboxplot(x[-which(x$MMRd_Type=="Delete"),], 
              x = z, y = y,
              #color = "Cancer_Type_NEW", 
              palette = "category20", 
              size = 0.5,
              add = c("mean","mean_sd"),
              add.params = list(color = "red")) +
              rotate_x_text(angle = 25) + 
              ylim(min_y, 2*max_y) + 
              ylab(ylabel) + 
              xlab("MMRd Type") +
              stat_compare_means(aes(label = ..p.signif..), comparisons = my_comparisons, label.y = 1.6*max_y*(1-0.02*(1:n)), vjust = 0.8, bracket.size = 0.1) + 
              stat_compare_means(label.y = 1.8*max_y) + 
              theme(axis.text = element_text(size = 10)) + 
              theme(axis.title = element_text(face = "bold")) + 
              theme(legend.text = element_text(size = 7)) + 
              theme(plot.margin = unit(c(0.5,0.5,0.5,0.5),"cm"))
 )
 return(get(paste0("p_",y)))    
  }else {
    max_y <- log10(max(x[,y]))
    min_y <- log10(min(x[,y]))
    assign(paste0("p_",y), ggboxplot(x[-which(x$MMRd_Type=="Delete"),], 
              x = z, y = y,
              #color = "Cancer_Type_NEW", 
              palette = "category20", 
              size = 0.5,
              add = c("mean","mean_sd"),
              add.params = list(color = "red")) + 
              #add = "jitter") + 
              rotate_x_text(angle = 25) + 
              ylim(min_y, 2*max_y) + 
              ylab(ylabel) + 
              xlab("MMRd Type") +
              stat_compare_means(aes(label = ..p.signif..), comparisons = my_comparisons, label.y = 1.6*max_y*(1-0.03*(1:n)), vjust = 0.8, bracket.size = 0.1) + 
              stat_compare_means(label.y = 2*max_y) + 
              scale_y_log10() +
              theme(axis.text = element_text(size = 10)) + 
              theme(axis.title = element_text(face = "bold")) + 
              theme(legend.text = element_text(size = 7)) + 
              theme(plot.margin = unit(c(0.5,0.5,0.5,0.5),"cm"))
 )
 return(get(paste0("p_",y)))
  }
  
}

# c(max_y*0.7, max_y*0.8, max_y*0.9)
metrics <- c("MSI_Score_1", "Impact_TMB_Score", "Fraction_Genome_Altered", "No_Missense", "No_Frameshift", "No_Splice","No_INSERT", "No_DELETION", "No_INFRAME","No_INDEL","No_INFRAME_Ins","No_INFRAME_Del","No_Nonsense","No_Frameshift_Ins","No_Frameshift_Del" )

df_2_pls_1 <- cbind(df_2[,! names(df_2) %in% c("Fraction_Genome_Altered","Cancer_Type_NEW","MMRd_Type","Row.names","Sample_ID","Tumor_Purity","Multiplex_or_not","MSH6_or_not","PMS2_or_not","Primary_Site","Metastatic_Site","Sample_Type","Race","Ethnicity")]+1, df_2[,c("Fraction_Genome_Altered","Cancer_Type_NEW","MMRd_Type","Tumor_Purity","Multiplex_or_not","MSH6_or_not","PMS2_or_not","Primary_Site","Metastatic_Site","Sample_Type","Race","Ethnicity")])

for(j in c("Multiplex_or_not","MSH6_or_not","PMS2_or_not")){
  for(i in metrics) {
  #compare means between MMRd groups, grouped by cancer types
  corr_tmp <- compare_means(as.formula(paste0(i," ~ ",j)), data = df_2[!df_2$MMRd_Type %in% c("?","Delete","Complex","NORMAL"),], group.by = "Cancer_Type_NEW")
  write.table(corr_tmp,
              file = file.path(getwd(),paste0(i,".corr.result.",j,".grouped.by.cancer.txt")),
              quote = FALSE,
              sep = "\t")

  # plot and test without grouping by cancer types
  assign(paste0("p_",i), plotCorr.MMRd(df_2_pls_1,i,j))
  print(get(paste0("p_",i)))
  
  ## Build a linear model for each genomic feature using predictors 1) MMRd type; 2) Cancer type; 3) Tumor purity
  sink(paste0("summary_lm_for_",i,"_w.",j,".Purity.txt"))
  model <- lm(reformulate(termlabels = c(j, 'Tumor_Purity'), response = i), df_2[df_2$Multiplex_or_not != "Delete",])
  print(summary(model))
  sink()
  sink(paste0("summary_annova_for_",i,"_w.",j,".Purity.txt"))
  print(anova(model))
  sink()
  sink(paste0("Summ.summary_lm_for_",i,"_w.",j,".Purity.txt"))
  print(summ(model, confint = TRUE, digits = 3))
  sink()
  p0 <- plot_summs(model, plot.distributions = T)
  print(p0 + ggtitle(paste0("Linear model: ", i, " ~ ", j, " + Tumor_Purity")) +
         scale_y_discrete(limits=levels(p0$data$term), labels=gsub(j,"",levels(p0$data$term))) + 
          theme(plot.title = element_text(size = 12, face = "bold")))
  
  sink(paste0("summary_lm_for_",i,"_w.",j,".txt"))
  model1 <- lm(reformulate(termlabels = c(j), response = i), df_2[df_2$Multiplex_or_not != "Delete",])
  print(summary(model1))
  sink()
  sink(paste0("summary_annova_for_",i,"_w.",j,".txt"))
  print(anova(model1))
  sink()
  sink(paste0("Summ.summary_lm_for_",i,"_w.",j,".txt"))
  print(summ(model1, confint = TRUE, digits = 3))
  sink()
  p1 <- plot_summs(model1, plot.distributions = T)
  print(p1 + ggtitle(paste0("Linear model: ", i, " ~ ", j)) +
          scale_y_discrete(limits=rev(levels(p1$data$term)), labels=gsub(j,"",rev(levels(p1$data$term)))) + 
          theme(plot.title = element_text(size = 12, face = "bold")))
  
  #Compare the two models w/ or w/o purity
  # export_summs(model, model1, error_format = "[{conf.low}, {conf.high}]",
  #              to.file = "xlsx",
  #              file.name = file.path(getwd(),paste0("model_",j,"_cmp_wo_purity_for_",i,".xlsx"))
  #              )
  # 
  
  # sink(paste0("summary_lm_for_",i,"_w.Purity.txt"))
  # model <- lm(reformulate(termlabels = c('Tumor_Purity'), response = i), df_1)
  #print(summ(model))
  }
}
## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Warning: Using `as.character()` on a quosure is deprecated as of rlang 0.3.0.
## Please use `as_label()` or `as_name()` instead.
## This warning is displayed once per session.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Loading required namespace: broom.mixed
## Registered S3 method overwritten by 'broom.mixed':
##   method      from 
##   tidy.gamlss broom
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.

## Warning: Removed 4 rows containing missing values (geom_segment).
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0.845098040014257, 1, 1.20411998265592, :
## cannot compute exact p-value with ties
## Warning in wilcox.test.default(c(0.845098040014257, 1, 1.20411998265592, :
## cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0.778151250383644, 0.698970004336019,
## 0.698970004336019, : cannot compute exact p-value with ties
## Warning in wilcox.test.default(c(0.778151250383644, 0.698970004336019,
## 0.698970004336019, : cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0.477121254719662, 0, 0,
## 0.477121254719662, : cannot compute exact p-value with ties
## Warning in wilcox.test.default(c(0.301029995663981, 0.301029995663981,
## 0.477121254719662, : cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0.845098040014257, 1, 1.20411998265592, :
## cannot compute exact p-value with ties
## Warning in wilcox.test.default(c(0.845098040014257, 1, 1.20411998265592, :
## cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0.301029995663981, :
## cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0.477121254719662, 0, 0,
## 0.301029995663981, : cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning in wilcox.test.default(c(0.778151250383644, 0.698970004336019,
## 0.698970004336019, : cannot compute exact p-value with ties
## Warning in wilcox.test.default(c(0.778151250383644, 0.698970004336019,
## 0.698970004336019, : cannot compute exact p-value with ties
## If some of the distribution curves are too short to see, consider rescaling
## your model coefficients or using the rescale.distributions = TRUE argument.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.

## Warning: Removed 2 rows containing missing values (geom_segment).
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.

## Warning: Removed 2 rows containing missing values (geom_segment).
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.ymin` is deprecated. Use `fun.min` instead.
## Warning: `fun.ymax` is deprecated. Use `fun.max` instead.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.

## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.